#pragma once

typedef struct
{
          int pos[4];
          complex entr[4];
} dirac_matr;

dirac_matr base_gamma[22];
dirac_matr Pplus, Pminus;

//Initialize a dirac matrix with outside entries
void init_dirac(dirac_matr *out,int pos0,double rea0,double ima0,int pos1,double rea1,double ima1,int pos2,double rea2,double ima2,int pos3,double rea3,double ima3)
{
          out->pos[0]=pos0;
          out->pos[1]=pos1;
          out->pos[2]=pos2;
          out->pos[3]=pos3;

          out->entr[0]=rea0 + I*ima0;
          out->entr[1]=rea1 + I*ima1;
          out->entr[2]=rea2 + I*ima2;
          out->entr[3]=rea3 + I*ima3;

}

void init_base_gamma(char *gamma)
{
          const double rad2=1./sqrt(2);
	  if (strcmp (gamma,"TWISTED") == 0)
	  {
		printf("TWISTED gamma basis\n");
          // TW convention
          init_dirac(&(base_gamma[0]),  0,1,0  , 1,1,0  , 2,1,0  , 3,1,0 );
          init_dirac(&(base_gamma[1]),  3,0,-1 , 2,0,-1 , 1,0,1  , 0,0,1 );
          init_dirac(&(base_gamma[2]),  3,-1,0 , 2,1,0  , 1,1,0  , 0,-1,0);
          init_dirac(&(base_gamma[3]),  2,0,-1 , 3,0,1  , 0,0,1  , 1,0,-1);
          init_dirac(&(base_gamma[4]),  2,-1,0 , 3,-1,0 , 0,-1,0 , 1,-1,0);
          init_dirac(&(base_gamma[5]),  0,1,0  , 1,1,0  , 2,-1,0 , 3,-1,0);
          init_dirac(&(base_gamma[6]),  1,0,1  , 0,0,1  , 3,0,1  ,  2,0,1);
          init_dirac(&(base_gamma[7]),  1,1,0  , 0,-1,0 , 3,1,0  ,  2,-1,0);
          init_dirac(&(base_gamma[8]),  0,0,1  , 1,0,-1 , 2,0,1  ,  3,0,-1);
          init_dirac(&(base_gamma[9]),  1,0,-1 , 0,0,-1 , 3,0,1  ,  2,0,1);
          init_dirac(&(base_gamma[10]), 1,-1,0 , 0,1,0  , 3,1,0  ,  2,-1,0);
          init_dirac(&(base_gamma[11]), 0,0,-1 , 1,0,1  , 2,0,1  ,  3,0,-1);
          init_dirac(&(base_gamma[12]), 3,0,1  , 2,0,1  , 1,0,1  ,  0,0,1);
          init_dirac(&(base_gamma[13]), 3,1,0  , 2,-1,0 , 1,1,0  ,  0,-1,0);
          init_dirac(&(base_gamma[14]), 2,0,1  , 3,0,-1 , 0,0,1  ,  1,0,-1);
          init_dirac(&(base_gamma[15]), 2,1,0  , 3,1,0  , 0,-1,0 ,  1,-1,0);

	  //sigma_dual Hay un signo - con las gammas chirales
          init_dirac(&(base_gamma[16]), 1,0,1 , 0,0,1 , 3,0,-1  ,  2,0,-1 );
          init_dirac(&(base_gamma[17]), 1,1,0 , 0,-1,0  , 3,-1,0  ,  2,1,0 );
          init_dirac(&(base_gamma[18]), 0,0,1 , 1,0,-1  , 2,0,-1  ,  3,0,1 );
          init_dirac(&(base_gamma[19]), 1,0,-1  , 0,0,-1  , 3,0,-1  ,  2,0,-1 );
          init_dirac(&(base_gamma[20]), 1,-1,0  , 0,1,0 , 3,-1,0  ,  2,1,0 );
          init_dirac(&(base_gamma[21]), 0,0,-1  , 1,0,1 , 2,0,-1  ,  3,0,1 );

            init_dirac(&Pplus ,          0,rad2,rad2  , 1,rad2,rad2  , 2,rad2,-rad2 , 3,rad2,-rad2);
            init_dirac(&Pminus,          0,rad2,-rad2 , 1,rad2,-rad2 , 2,rad2,rad2  , 3,rad2,rad2 );
	}

          if (strcmp (gamma,"APE") == 0)

          {
		printf("APE gamma basis\n");

          init_dirac(&(base_gamma[0]),  0,1,0  , 1,1,0  , 2,1,0  , 3,1,0 );
          init_dirac(&(base_gamma[1]),  3,0,1  , 2,0,1  , 1,0,-1 , 0,0,-1);
          init_dirac(&(base_gamma[2]),  3,1,0  , 2,-1,0  , 1,-1,0  , 0,1,0);
          init_dirac(&(base_gamma[3]),  2,0,1 , 3,0,-1  , 0,0,-1  , 1,0,1);
          init_dirac(&(base_gamma[4]),  0,1,0 , 1,1,0 , 2,-1,0 , 3,-1,0);
          init_dirac(&(base_gamma[5]),  2,1,0  , 3,1,0  , 0,1,0 , 1,1,0);

          init_dirac(&(base_gamma[6]),  1,0,1  ,  0,0,1 , 3,0,1 , 2,0,1  );
          init_dirac(&(base_gamma[7]),  1,1,0  , 0,-1,0 , 3,1,0 , 2,-1,0 );
          init_dirac(&(base_gamma[8]),  0,0,1  , 1,0,-1 , 2,0,1 , 3,0,-1);
          init_dirac(&(base_gamma[9]),  3,0,1  , 2,0,1  , 1,0,1 , 0,0,1);
          init_dirac(&(base_gamma[10]), 3,1,0  , 2,-1,0 , 1,1,0 , 0,-1,0);
          init_dirac(&(base_gamma[11]), 2,0,1  , 3,0,-1 , 0,0,1 , 1,0,-1 );
          init_dirac(&(base_gamma[12]), 1,0,1  , 0,0,1  , 3,0,-1, 2,0,-1 );
          init_dirac(&(base_gamma[13]), 1,1,0  , 0,-1,0 , 3,-1,0, 2,1,0 );
          init_dirac(&(base_gamma[14]), 0,0,1  , 1,0,-1 , 2,0,-1, 3,0,1);
          init_dirac(&(base_gamma[15]), 2,1,0  , 3,1,0  , 0,-1,0, 1,-1,0 );

	  //sigma_dual

          init_dirac(&(base_gamma[16]), 3,0,1  , 2,0,1  , 1,0,1 , 0,0,1 );
          init_dirac(&(base_gamma[17]), 3,1,0  , 2,-1,0 , 1,1,0 , 0,-1,0);
          init_dirac(&(base_gamma[18]), 2,0,1  , 3,0,-1 , 0,0,1 , 1,0,-1  );
          init_dirac(&(base_gamma[19]), 1,0,1  ,  0,0,1 , 3,0,1 , 2,0,1  );
          init_dirac(&(base_gamma[20]), 1,1,0  , 0,-1,0 , 3,1,0 , 2,-1,0  );
          init_dirac(&(base_gamma[21]), 0,0,1  , 1,0,-1 , 2,0,1 , 3,0,-1 );


	}


}

int sigma_dual(int sigma)
{
	int g5sig;
	 if (sigma==6) g5sig=9;
         else if (sigma==7) g5sig=10;
         else if (sigma==8) g5sig=11;
         else if (sigma==9) g5sig=6;
         else if (sigma==10) g5sig=7;
         else if (sigma==11) g5sig=8;
	 else g5sig=0;
	 return g5sig;

}

